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1 Introduction 

Let A, B, C £ R" xn be given matrices, the generic Quadratic Eigenvalue Comple¬ 
mentarity Problem (QEiCP) widely studied in recent papers, e.g., see [6,14,38], is 
that of finding (A, ijelx R" such that 



K 9 a:_L(A 2 A + A B + C)x £ K*, 


( 1 . 1 ) 


where e= (1,1,..., 1) T £ R", K is a closed convex cone in R", and K* refers to 
the dual cone of K, which is defined by 


k* := {y ^ I (y,x) > 0, V x £ K} . 

Without loss of generality, the linear constraint e T x = 1 in (1.1) plays an impor¬ 


tant role in preventing the x component of a solution to vanish. Notice that the 
leading matrix A could be singular, and in particular, the QEiCP immediately re¬ 
duces to the classical Eigenvalue Complementarity Problem (EiCP) when A = 0. 
Clearly, QEiCP is an interesting generalization of the classical EiCP with embrac¬ 
ing an extra quadratic term on A. When K is taken as the whole space R", then 
(1.1) becomes the well studied unconstrained quadratic eigenvalue problem, which 
frequently arises in areas such as the electric power systems [28] , the dynamic anal¬ 
ysis of acoustic systems [3] and linear stability of flows in fluid mechanics [27], to 
name a few. We refer the reader to [41] for a survey on the unconstrained version. 
If the matrices A, B, C are all symmetric, QEiCP and EiCP are called symmetric, 
respectively. 

Since the seminal work on EiCP [10] devoted to the study of static equilibrium 
states of mechanical systems with unilateral friction, both EiCP and QEiCP have 
been well discussed from theoretical and numerical perspective in the literature, 
e.g., see [1,11,14,15,16,22,23,24,38], where these papers only focus on matrix 
cases. For instance, in order to study the sufficient conditions for the existence of 
solutions of QEiCP, the so-called co-regularity and co-hiperbolicity properties were 
introduced in [38]. Usually, the co-regularity on matrix A £ R" xn , i.e., x T Ax i=- 0 
for any x £ K , means that A or —A is strictly A-positive. We say that ( A , B, C) € 
M n := R" xn x R nx ” x R nxn satisfies co-hiperbolicity property, if 


(x 1 Bx) 2 > 4(;r T Ax)(x T Cx), M x £ K. 


However, checking whether a given matrices triplet ( A,B,C ) £ M n satisfies co- 
hyperbolicity or not is co-NP-complete, which is essentially the verification prob¬ 
lem of copositiveness (see Definition 2.1) of a related fourth order n-dimensional 
tensor [32,39]. 
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In recent decade, tensor, which is a natural extension of the concept of ma¬ 
trices, is on the timely topic of high-dimensional data representation in terms of 
theoretical analysis and algorithmic design because of its widespread applications 
in engineering. A tensor, namely, is a multidimensional array, whose order is the 
number of dimensions. Let m and n be positive integers. We call A = {a i±im ), 
where £ R for 1 < < n, a real m-th order n-dimensional real 

square tensor. For the sake of convenience, we denote by T m ,n the space of m-th 
order n-dimensional real square tensors. Furthermore, a tensor A £ Tm,n is called 
symmetric if its entries are invariant under any permutation of its indices. For a 
vector x = (xi ,..., x n ) T £ C n and a tensor A = £ T m ,n, Ax 171-1 is an 

n-dimensional vector with its i-th component defined by 

n 

(Ax )j ^ ^ 1, 2,..., n, 

and Ax m is the value at x of a homogeneous polynomial, defined by 

n 

Ax — ^ 2 . ..i rn Xi 1 Xi 2 • • • Xi m . 

Although tensor-related problems have been received considerable attention 
many years ago, the history of research on eigenvalues (eigenvectors) of a square 
tensor can be traced back to the pioneer works independently introduced by Qi 
[31] and Lim [25]. Comparatively speaking, the developments of eigenvalue-related 
problems for tensors are still in their infancy. For given tenors A,B £ Tm,n, we say 
that (A, B) is an identical singular pair, if 

|x £ C n \{0} | Ax’ 71-1 = 0, Bx 771 ” 1 = oj / 0. 

When assuming that ( A , B) is not an identical singular pair, we say (A, x) £ C x 
(C"\{0}) is an eigenvalue-eigenvector pair of (A,B), if the following n-system of 
equations 

(A-XB)x m ~ 1 =0 (1.2) 

possesses a nonzero solution x. This unified definition of eigenvalue-eigenvector 
pair for tensors was introduced by Chang et al. [8]. In recent years, it is well docu¬ 
mented in the literature that tensors and eigenvalues/eigenvectors of tensors have 
fruitful applications in various fields such as magnetic resonance imaging [4,34], 
higher-order Markov chains [29] and best-rank one approximation in data analy¬ 
sis [33], whereby many nice properties such as the Perron-Frobenius theorem for 
eigenvalues/eigenvectors of nonnegative square tensors have been well established, 
see, e.g., [7,43]. All these encourage us to consider tensor eigenvalue complemen¬ 
tarity problems. However, to the best of our knowledge, the most recent paper [26] 
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is the first work devoted to the Tensor Generalized Eigenvalue Complementarity 
Problem (TGEiCP), whereas leaving higher-degree cases a big gap. Therefore, the 
main objective of this paper is to £11 out this gap. 

In this paper, we consider the Tensor Higher-Degree Eigenvalue Complemen¬ 
tarity Problem (THDEiCP), which goes beyond the framework of QEiCP and 
further generalizes TGEiCP. Mathematically, the THDEiCP can be characterized 
as finding a scalar A € R and a vector x £ R”\{0} such that 

K B x A. {\ m A +XB+ C)x m ~ 1 £ A*, (1.3) 

where A = (ai 1 i 2 ...i m ), B = andC = (c ili2 ... irn ) £ T m ,n- Correspondingly, 

the scalar A and the nonzero vector x satisfying system (1.3) are respectively 
called an m-degree A'-eigenvalue of the tensors triplet Q := (A,B,C) £ AFm,n := 
T m ,n x Tm,n x T m ,n and an associated A-eigenvector. Alternatively, (A,*) is also 
called an m-degree A'-eigenpair of Q. Throughout, the set of all m-degree K- 
eigenvalues of Q is called the m-degree A-spectrum of Q, i.e., 

a(Q, K) := |a £ R | 3 x £ R"\{0}, A 9 1 1 (X m A + XB + C)x m_1 £ A*j . (1.4) 

If A : = R” , the m-degree A'-eigenvalue/eigenvector of the tensors triplet Q is called 
the m-degree Pareto-eigenvalue/eigenvector of Q, and the m-degree A-spectrum 
of Q is called the m-degree Pareto-spectrum of Q. If * £ int(A) (resp. x £ R"+), 
then A is called a strict m-degree A'-eigenvalue (resp. Pareto-eigenvalue) of Q. 

It is clear from (1.3) that THDEiCP covers TGEiCP and QEiCP as the special 
cases. More concretely, by taking A = 0, model (1.3) immediately reduces to 
TGEiCP studied in [26]. When we set m = 2, THDEiCP clearly corresponds to the 
QEiCP (1.1). Like [37], on the other hand, we can further establish the connection 
between THDEiCP and a class of differential inclusions with nonconvex processes 
F defined by 

Gr(r) := {(x,y) £ E" x 1" | A 9 x±{Ay m ~ 1 - Bx m ~ x ) £ A*}. 

Accordingly, for the differential inclusions defined by u(t) £ A(u(t)), as noticed 
already by Rockafellar [35], the change of variables u(t) = exp(At)* with A > 0 
leads to the equivalent system Xx £ F(x). Therefore, if the pair (A,*) satisfies Xx £ 
F(x), then the trajectory t 1-4 exp(At)* is a solution to the considered differential 
inclusions. Moreover, if the trajectory constructed above is nonconstant, then x 
must be a nonzero vector; this requires (A, x) to be a solution of THDEiCP with 
C = 0 because of A > 0. 

The paper is divided into six sections. As far as we know, it might be the first 
work on THDEiCP, we thus do not know whether the topological properties of 
QEiCP still hold for the newly introduced model. In Section 2, we first study, 
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as briefly as possible, some topological properties such as closedness , boundedness , 
and upper-semicontinuity of the m-degree Tv-spectrum given by (1.4), in addition 
to estimating upper bounds on the number of eigenvalues of THDEiCP. With the 
preparations on these topological properties, ones may be further concerned with 
that how to solve the model under consideration. In Section 3, we reformulate the 
special case of THDEiCP with symmetric A and B as a weakly coupled polyno¬ 
mial optimization problem for the case where K := R” , which would potentially 
facilitate the algorithmic design. From a theoretical perspective, in Section 4, we 
establish the results concerning existence of solutions of THDEiCP without sym¬ 
metry assumptions on A and B. Based upon the augmented Lagrangian method, 
in Section 5, we propose an implementable splitting algorithm to solve the result¬ 
ing polynomial optimization reformulation of the symmetric THDEiCP and report 
some preliminary results. Finally, we give some concluding remarks in Section 6. 

Notation. Let R™ denote the real Euclidean space of column vectors of length n, 
which is equipped with the standard inner product (y, x) = y T x and the associated 
norm. The superscript ‘T’ indicates transposition and the symbol ‘_L’ represents 
orthogonality. Denote R" = {x £ R" | x > 0} and R" + = (x € R n | x > 0}. For 
a vector x £ R" and a subset J of the index set [n] := {1,2,... , n}, we use the 
notation xj for the | J| dimensional sub-vector of x, which is obtained by deleting 
all components i qL J, where the symbol |J| denotes the cardinality of J. For a 
vector x £ R" and an integer r > 0, denote = (x\,X 2 , ...,xh) T , and denote 

by diag(*) the n x n diagonal matrix containing xt in its diagonal. Moreover, for 
A £ 7m,n, we denote the principal sub-tensor of A by Aj, which is obtained by 
homogeneous polynomial Ax m for all x = ( x\,X 2 , ■ ■ ■ ,x n ) T with Xi = 0 for [n]\J. 
So, „4j £ T m jjp Denote by e £ R" with all entries being 1, i.e., e = (1, 1,..., 1) T . 
Denote by X = (<5q...i m ) the unit tensor in Tm,n, where is the Kronecker 

if u = . . . = imi 
otherwise. 



2 Some basic properties of iVT-spectrum 

In this section, we summarize some basic definitions and study some basic topo¬ 
logical properties of m-degree cone spectrum, which will be used in subsequent 
sections. 

We first give the concept of cone positive square tensors, which is a generalized 
concept of copositive square tensor introduced in [32] and studied in [39]. 
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Definition 2.1 Let K be a closed convex cone in R™ and Q £ T m ,n- We say that Q 
is a (resp. strictly) A'-positive tensor, if Qx m > 0 (resp. > 0) for any x £ A\{0}. If 
K = R”, the (strictly) A'-positive tensor Q is said the (strictly) copositive tensor. 

It is obvious from the notation of T m ,n that T m ,n is a linear space. The distance 
between two elements Qi = ( Ai,Bt,Ci ) G -Fm,n (i = 1,2) is measured by means of 
the expression 

l!Qi - Q 2 \\f = {ll-Ai -A 2 \\ 2 f + ||Bi - B 2 \\f + \\Ci -C 2 |||} 3 , 

where 

\\A\\F = / y -A.= (fli 1 ...i Tn ) G 7m,n- 

Denote by ^f(R n ) the set of nonzero closed convex cones in R”, which is asso¬ 
ciated with the natural metric defined by 

S(Ki, K 2 ) := sup |dist(z, A'i) — dist(z, K 2 )\, 

Ml<i 

where dist(«, A') := inf. ug ^-||« — u|| stands for the distance from 2 to A'. An equiv¬ 
alent way of defining S is 

5(Ki,K 2 ) = haus(AT n ® n , K 2 n ®„), 

where ® n is the closed unit ball in R n , and 

haus(^i,^ 2 ) :=max-j sup dist(a, %), sup dist(z,'^i) 

stands for the Hausdorff distance between the compact sets r £\, r fi 2 C R™ (see [2, 
pp. 85-86]). General information on the metric 5 can be consulted in the book by 
Rockafellar and Wets [36] . According to [42] , the operation K i-)- K* is an isometry 
on the space (if (R”), (5), that is to say, 

5{Kf,K 2 ) = S(K 1 ,K 2 ), for all A'i, A 2 £ V(R n ). 

The basic topological properties of the mapping a : T rn .n x ^(R") —» 2 R , defined 
in (1.4), are listed in the next proposition. This proposition is a tensor version of 
generalizing the results presented in [38]. As far as semicontinuity concepts are 
concerned, we use the following terminology (cf. Section 6.2 in [2]): 

Definition 2.2 Let W and Y be two topological spaces. The mapping <P : W —> 2 y 
is said to be upper-semicontinuous (resp. lower-semicontinuous) if the set 

{w£W \ A(w) C U} (resp. {w £ W \ &(w) D U ^ 0}) 

is open, whenever U C Y is open. 
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Definition 2.3 Let Q = ( A,B,C) £ F rn ,n and K £ if(R n ). We say that Q is 
A-regular if the leading tensor A satisfies 

Ax m ^ 0, Vx£ 7A{0}. 

It is obvious that, if Q is A'-regular, then the leading tensor A in Q or —A is 
A'-positive. 

Proposition 2.1 The following three statements are true: 

(i) . The set £ := {(Q,K, A) £ Tm,n x ^(K n ) x R | A £ <x{Q,K)} is closed in the 
product space T m ,n x ^(R 11 ) x R. In particular, for any ( Q,K ) £ T m ,n x ^(R™), 
a(Q,K) is a closed subset of R; 

(ii) . Let ( Q,K ) £ J- m ,n x ^(R™). If Q is K-regular, then the mapping a : J- m ,n x 

^(R”) —>■ 2 r is locally bounded at ( Q,K ), i.e., U(Q K)eN bounded for 

some neighborhood AT of (Q,K). 

(iii) . Let (Q,K) £ T m , n x ^(R"). If Q is K-regular, then a is upper-semicontinuous 
at (Q, K). 

Proof (i). The closedness of £ amounts to saying that 

j 

Since £ o-(Qv, AT), there exists a vector xj, £ R n \{0} such that 

AT 3x„l (A™Av + \vB v + tT)®™” 1 £ A'*. (2.1) 

Let 5T = x v l\x v |. From the homogeneity of (2.1) on x , we know that 

AT (A™A, + A„Z3„ + C v )x™~ x £ AT*. (2.2) 

It is clear that ||x„|| = 1 for every v = 1,2,.... Without loss of generality, we 
assume that x v —> x. It is obvious that ||x|| = 1, which means x £ R™\{0}. Since 
S(Ki, A'J) = <5(AT, AT) for any AT, AT £ ^(R n ), by passing to the limit in (2.2), 
one knows 

Kbx 1 (A m A + A B + C)x m ~ 1 £ K*, 

which implies A £ cr(Q, K) due to x £ R”\{0}. We proved the first part (i) of this 
proposition. 

(ii). Suppose that the map a is not locally bounded at ( Q, K ). Then there 
exists a sequence {Q V ,K V ,X V } satisfying 

\\Qv-Q.\\f^Q, 5(AT,A')^0, and |A„| -»■ + 00 , 
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such that A„ £ cr(Qv,K v ) for any v = 1,2,.... Consequently, there exist vectors 
x v £ K v with ||a;i/|| = 1, such that 

K v 3x„± (A ™A V + A V B V + £ K* v . (2.3) 


By (2.3), we have 


(A™ Au + XisBu + C v )x™ = 0, 


which implies 



B v C v \ 
a ™- 1 a™; 


= 0 . 


We assume, without loss of generality, that x v —> x. It is obvious that ||x|| = 1, 
which means x £ R™\{0}. By passing to the limit in the above expression, it holds 
that Ax m = 0. It contradicts the If-regularity of Q by the truth * £ I\\{0}. 

(iii). Suppose that a is not upper-semicontinuous at ( Q, K ), then we could find 
an open set U C R and a sequence {(Qi/, A'„)} satisfying (Qj,,/W) —> ( Q, K ), such 
that 


a(Q,K)cU but a(Q v ,K v )n{R\U) ^ 0, for any v = 1, 2, .... 

Now, for each v, pick up \ v £ a(Q^,K u ) n (R\f/). It follows from (ii) that the 
sequence {A^} admits a converging subsequence. By (i), the corresponding limit 
must be in a(Q,K ) n (W\U), which together with a(Q,K) C U leads to a contra¬ 
diction. □ 

From the first two parts (i) and (ii) of Proposition 2.1, we have the following 
corollary. 

Corollary 2.1 Let (Q, I\) £ J- m ,n x If Q is K-regular, then a(Q,K) is com¬ 

pact. 

Below, we present a preliminary estimation on the numbers of m-degree Pareto- 
eigenvalues. We first present the following proposition which fully characterizes the 
m-degree Pareto-spectrum of THDEiCP. 

Proposition 2.2 Let Q = ( A,B,C ) £ Tm,n- A rea l number A is an m-degree Pareto- 
eigenvalue of Q, if and only if there exists a nonempty subset J C [n] and a vector 
w £ such that 

(\ m Aj + A Bj + Cj)w m_1 = 0 (2.4) 

£ (A m + A6ji 2 ...j m + ) w i 2 ---w im >0, V i £ [n]\J. 

*2 J 


and 
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In such a case , the vector x £ R™ defined by 

) Wi, i £ J, 

X' — \ 

I 0, * £ [n]\J 

is a Pareto-eigenvector of Q, associated to the m-degree Pareto-eigenvalue X. 

Proof It can be proved in a similar way that used in [40]. □ 

It is well known that, on the left-hand side of (2.4), (A m Aj-|-A£>j-|-Cj)w m ~ 1 is 
indeed a set of | J| homogeneous polynomials in | J| variables, denoted by {P^(w) | 1 < 
i < | J|}, of degree (m — 1). In the complex field, in order to study the solution 
set of a system of | J| homogeneous polynomials (Pi,... ,P|j|), in | J| variables, the 
concept of the resultant Res(Pi,... ,P|j|) is well defined and introduced, we re¬ 
fer to [9] for details. Applying to our current problem, Res(Pi,..., P| j|) has the 
following properties. 

Proposition 2.3 We have the following results: 

(i) . Res(Pi,..., P\j\) = 0, if and only if there exists (X, x) £Cx(C^\{0}) satisfying 
the relation (2.4). 

(ii) . The degree of X in Res(Pi,..., P| u) is at most m| J|(m — 1)^ _1 . 

By Proposition 2.2, if A is an m-degree Pareto-eigenvalue of Q = ( A,B,C ) £ 
Pm,n, then there exists a nonempty subset J C [n] such that A is a strict m-degree 
Pareto-eigenvalue of Qj = ( Aj,Bj,Cj) £ We now state and prove one of 

main results in this section. 


Theorem 2.1 Let Q = ( A,B,C ) £ T m ,n- Assume that Q is Mff-regular. Then, there 
are at most T m ,n '■= nm n m-degree Pareto-eigenvalues of Q. 

Proof It is obvious that, for every k = 0,1,..., n — 1, there are ( n " fc ) corresponding 
principal sub-tensors triplet of order m dimension ( n—k ). Moreover, by Proposition 
2.3, we know that every principal sub-tensors triplet of order m dimension (n — k ) 
can have at most m(n — fc)(m — strict m-degree Pareto-eigenvalues. By 

Proposition 2.2, in this way one obtains the upper bound 

n— 1 / 

Tm.n — / I 

fe =oV 

We obtain the desired result. □ 


n 

n—k 


m(n — k)(m — 1) 


n—k —1 


Now we extend the above result to the more general case where A is a poly¬ 
hedral convex cone. A closed convex cone K in R" is said to be finitely generated 
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if there is a linearly independent collection {r?i, 772 , • • ■, Vp} of vectors in R" such 
that 

I\ = cone{??i,r?2, ...,r) p } = | a = (01, 02, ■ ■ • ,a P ) T £ R+j • ( 2 . 5 ) 

It is clear that K = {H T a \ a £ R(j_}, where H = [771 , 772 , ■ ■ ■ ,Vp] ■ Moreover, it is 
easy to see that the dual cone K* of K is equivalent to {z £ R n | Hz > 0}. 

Theorem 2.2 Let Q = (A,B,C) £ T m ,n- Let K be represented by (2.5). Assume that 
Q is K-regular. Then, there are at most T m ,p '■= pmP m-degree K-eigenvalues of Q. 

Proof We first prove that problem (1.3) with K defined by (2.5) is equivalent to 
finding a vector d £ R p \{0} and A £ R such that 

a > 0, (\ m T> + \Q + 5)a m_1 > 0, (d, (A m H + A0 + 5)a m_1 ) = 0, (2.6) 


where V, Q and S are three m-th order p-dimensional tensors, whose elements are 
denoted by 

n 

a jij2-- j m hi 1 j 1 hi 2 j 2 ■ ■ ■ hi m j m , 


9iii 2 . 


E 


u Ji 32 


,hi 1 j 1 hi 2 j 2 • 




■ dm =1 


and 

n 

Si^i 2 ...i rn = ^ • ■■ hi m j m 

ji = l 

for 1 < *i,* 2 , • ■ ■, im. < P, respectively. 

Let ( X,x ) £ R x (R"\{0}) be an m-degree Iv-eigenpair of Q. Since x £ X\{0}, 
by the definition of K, there exists a nonzero vector a £ R p such that x = H T a. 
Consequently, from the fact that (A m ^l + A B + C)i m ~ 1 £ K* and the expression 
of K* , it holds that H(X m A + A B + C)x m ~~ x > 0, which implies 

H(X m A + \B + C){H T a) m ~ 1 >0. (2.7) 


By the definitions of V, Q and S, we know that (2.7) can be equivalently written 
as 

(X m D + A g + S) d m_1 > 0. 

Moreover, it is easy to verify that (d, (\ m V + A Q + 5)d m_1 ) = 0. Conversely, if 
(A,d) £ R x (R p \{0}) satisfies (2.6), then we can prove that (A, x) with x = H T a 
satisfies (1.3) in a similar way. 

Since Q is A'-regular, it is easy to verify that ( [V,Q,S ) is R^-regular. Conse¬ 
quently, by applying Theorem 2.1 to problem (2.6), we know that Q has at most 
T m ,p = pm p m-degree Tv-eigenvalues. □ 
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The above theorem shows that a(Q,K) has finitely many elements in case 
where K is a polyhedral convex cone. However, the situation can be even worse 
in the nonpolyhedral case. For instance, Iusem and Seeger [21] successfully con¬ 
structed a symmetric matrix C (i.e., Q = (0,1,— C) G Ti.n) and a nonpolyhedral 
convex cone K such that a(Q,K) behaves like the Cantor ternary set, i.e., it is 
uncountable and totally disconnected. 


3 Optimization formulation of THDEiCP 

In this section, for the purpose of finding solutions of THDEiCP, we introduce an 
optimization reformulation, which paves the way of designing algorithms. Here, 
we only consider the case where A and B are two symmetric tensors, C := —Z, and 
K := R". 

We consider the following homogeneous polynomial optimization problem. 


max <po(w, v) rn(m — 1) m 1 n T iJ m ^ — Bu m 
s.t. Au m +v r vi m -^ = 1, (3.1) 

u > 0, v > 0. 

Let 4>o(u,v) := Au m + v T v^ m ~ 1 \ and let <po be defined by (3.1). We derive that 

Vu^po(u,v) = m(m — 1) v* diag(v)iJ m ~ 2 l — mBu m ~ 1 , (3.2a) 

X? v <po(u,v) = m(m — 1) ™ _1 rJ rra_1 ], (3.2b) 

< 

Vtt<)>o(M, v) = viAu 771 - 1 , (3.2c) 

, V v (j>o(u, v) = rnr^ m_1 l (3.2d) 


Now, we state the relationship between (1.3) and (3.1) as follows. 

Theorem 3.1 Let Q = (A, B, —Z) G Assume that A and B are both symmetric. 

Let (u,v) withu 0 be a stationary point of problem (3.1). Then (A, u) is an m-degree 
Pareto-eigenpair of Q, where A = (ipo(u, v )) m ~ 1 . 

Proof Since (u,v) is a stationary point of problem (3.1), it follows from (3.2) that 


there exist a,/3 G R™ and 7 G R, such that 

m.Bu nl — m(m — 1) "‘diag^)^™” 2 ! = a + ^mAu n 1 , (3.3a) 

—m(m— l)™ _1 «^ n_1 ^ =(0 + 7 (3.3b) 

* a > 0, u > 0, a T u = 0, (3.3c) 

/3 > 0, v > 0, f) T v = 0, (3.3d) 

. Au m + v T v [m ~ 1] = 1. (3.3e) 
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Rearranging terms of (3.3b) yields 

— (m — 1 ) = 0/m. (3-4) 

We claim that 0 = 0. Otherwise, if 0 j=- 0, then there exists an index £ [n] 
such that 0i o > 0, which implies Vi 0 = 0 from (3.3d), and hence, it holds that 
— (m — l)m - 1 M ™ -1 = 0i o /m > 0. It is a contradiction. Therefore, it follows from 
(3.4) that 

— (m — 1) ^ (3-5) 

Moreover, it is clear from the truth u 0= 0 and (3.5) that 7 < 0 and 

v = (— 7 )” m -i (m — 1 )~™u. (3-6) 

By invoking (3.3a) and (3.6), we have 

Bu m ~ x - (- 7 )-str u [m ~ 1] = - > 0, 

m 

which implies 

(- 7)^1 Au m ~ x + ^)^Bu m ~ x - u [m ~ 1] > 0. (3.7) 

Moreover, using (3.6), (3.3a), and (3.3c), it is not difficult to verify that 

<«, = 0. (3.8) 

On the other hand, it follows from (3.3a) and (3.3c) that 

mBu m — m'/Au n — m(m — 1)"*« T «^ - ^ = 0, 

which implies 

—m<po(u, v) = m~/Au rn — m(m — 1 ) ^~ 1 v T v}- m ^ 

- a —m - —T -\m— ll 

= mjAu + mjv v l 1 

= my, (3.9) 

where the second equality comes from (3.5), and the last equality is due to (3.3e). 
Hence, we conclude from (3.9) that ipo(u,v) = —7 > 0, and both (3.7) and (3.8) 
mean that u is an eigenvector of (3.1) associated to the eigenvalue A. □ 

Theorem 3.2 Let Q = (A, B, —X) £ Assume that A and B are both symmetric 

and A is copositive. Let (A, S) be an m-degree Pareto-eigenpair of Q. If X > 0, then 
(u,v) is a stationary point of (3.1), where 

(u,v) = - ±{x,y) 

(. Ax m + y T y I™- 1 !) m 

with y = (m — 1) _ ™ (A) _ 1 i. 


(3.10) 
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Proof Since y £ R+\{0} and A is copositive, we have Ax m + > 0. More¬ 

over, it is easy to check that (u,v) given in (3.10) is a feasible solution of (3.1). 
Take 

a = ----^ (A m Ax m ~ 1 + A Bx m ~ x - , 

\(Ax m + J/ T y[ m_1 l ) m ' ' 

j3 = 0 and 7 = —A" 1-1 . It is obvious that d > 0, since (A, x) satisfies (1.3) and 
A > 0. Moreover, it is not difficult to see that a T u = 0 and j3 T v = 0. Finally, with 
the definition of (u,v) given in (3.10), we can verify that (3.3a) and (3.3c) hold, 
which means that (u,v) is a stationary point of (3.1). □ 

Denote w := (« T , v T ) T and = Wi for i = 1,2,, 2 n. When A is strictly 

copositive, the feasible set of (3.1) is compact. Hence, the globally optimal value 
of (3.1), denoted by ^™ ax , exists. Denote by w the corresponding globally optimal 
solution, and denote 

d = ((ei(w)) T >-f(w'j<=( t5 )) T ) T with t = ^ 

iG/ (w) 

where I(w) = {i £ [2 n] 1^ = 0} and I c (w) = [n]\J(to). From the homogeneity of <j>o 
and the fact that <j>o(w) = 1 , it is easy to verify that w T = rrup o(w) = m^ 0 , 

which implies that V</>o(w) 0 and hence it is linearly independent. Moreover, it 

is not difficult to know that 

d T V<f>o(w) = m I (Au n ~ 1 ) i — t I = 0 
\ieUw) ) 

and d T V0i(wi) = 1 > 0 for every i £ I(w). This means that the Mangasarian- 
Fromovitz constraint qualification (MFCQ) holds at w. Therefore, we know that 
w is a stationary point of (3.1). Moreover, we claim that ii^O. In fact, by taking 

( _ m \ 1/m 

j e with a = i m =% a*!...j m , we can see that 

wt = ( ut,vt .) is a feasible solution of (3.1), and the corresponding objective value 
is 

Po(u t , vt) = t m ^nm(m - 1)™ _1 - bt^j , 

wher eb = J2i 1 ,...,i m = 1 b i 1 -im- Hence, we know that <po(ut,vt) > 0 for t > 0 enough 
small, which implies that <po(w, v) > 0 due to the fact that (u,v) is an optimal 
solution of problem (3.1). Consequently, it holds that u^0. Moreover, by Theorem 
3.1, we know that (A ,u) with A = (ipo(u, v)) m ~ 1 is a solution of (1.3), which implies 
that (1.3) has at least a positive m-degree Pareto-eigenvalue. Therefore, one has 
Vo ' ax < Amax 1 , where 

Amax = max {A £ R | 3 x £ R”, (A,*) is an m-degree Pareto-eigenpair of Q} . 
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Theorem 3.3 Let Q = (A,B, —I) £ Assume that A and B are both symmetric 

and A is strictly copositive. Then, we have 

\ m— 1 _ max 

''max — TO 


Proof Let (A, x) with A > 0 be an m-degree Pareto-eigenpair of Q. By the homo¬ 
geneity of the complementarity system (1.3) with respect to x, without loss of 
generality, we assume that x satisfies x = 1. As a consequence, it immediately 
follows from the strict copositiveness of A that ^lr m > 0. Denote 


V = 


(m 



and (m, v) 


_ (x,y) _ 

( Ax m + y^yt m — 1 1 ) 


(3.11) 


It is trivial that (u, v) £ R" x R™ , 


Au m = 


Ax m -I- y T y\. m ~ x \ 


Ax r ‘ 


and v v 


-T-[m— 1] 

'll 1 J — 


Ax m + yTyl" 1 - 1 ! 


-T -\m— 1] 

y y [ , 


which implies that Au m + v T = 1 holds. Hence, ( u,v) is a feasible solution 
of (3.1) and ipo(u,v) < v J o iax - 

On the other hand, since (A, x) is an m.-degree Pareto-eigenpair of Q , we know 
that X m Ax m + A Bx m — x 1 = 0. Substituting (u, v) into ipo{u,v) yields 


ipo(u, v) = m(m — 1) m 1 u T u^ m ^ — Bu m 

_ mx T x^ m -^ - (m - 1)A Bx m - ^ 
(?n — 1)A m Ax m + 

== A m_1 , 


where the second equality comes from (3.11). Therefore, it holds that A m_1 < 
(p™ 3,51 , which implies that AJO^x 1 < (p™ 351 . To sum up, we obtain the desired result 
of this theorem. □ 


4 Existence of solutions for THDEiCP 

In this section, we study more general results on the existence of solutions of 
THDEiCP with C := ±1 and K = R”, but without symmetry assumptions on A 
and B. We first present the existence result of symmetric tensors. 

Theorem 4.1 Let Q = (A, B, —T) £ Tm,n- Assume that A and B are both symmetric 
and A is strictly copositive. Then Q has at least an m-degree Pareto-eigenpair. 

Proof Consider the homogeneous polynomial optimization problem (3.1). Since 
A is strictly copositive, it is easy to see that the feasible set of (3.1) is compact. 
Consequently, from the continuity of the objective function (po in (3.1), its globally 
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optimal solution, denoted by (u,v), exists. As arguments above, the constraint 
qualification MFCQ holds at (u,v). Hence, (u,v) is a stationary point of (3.1). 
Moreover, we know that ipo(u,v) > 0 from the arguments above, which implies 
that u 0. Therefore, the assertion of Theorem 3.1 shows that (A, u) is an m- 
degree Pareto-eigenpair of Q, where A = (ipo(u, v)) . □ 

The above theorem is a fundamental result for THDEiCP. However, many 
real-world problems often violate the symmetry condition. In other words, the 
symmetry assumptions on A and B are relatively stronger. Indeed, a general ex¬ 
istence theorem of solutions of QEiCP have been well established in [38], which 
states that, if (A, B, C) satisfies co-hyperbolicity properties and the leading matrix A 
is co-regular, then the considered QEiCP has at least a solution. As a generalization 
of QEiCP, we are naturally concerned with whether such a similar result of QE¬ 
iCP also holds for tensors. Hereafter, we study such a more general result without 
assuming the symmetry of A and B , in addition to presenting some checkable 
conditions on Q = ( A,B,C ) instead of the co-hyperbolicity. 

Theorem 4.2 Let Q = (A,B,—I) £ F m ,n. Assume that A is strictly copositive. If A 
and B satisfy the following condition 

{au...i + 1 - m)(m - 1)™ _1 - bu...i > 0, V i £ [n], (4.1) 

then Q has at least an m-degree Pareto-eigenpair. 

Proof We first denote two sets by 

S = {{x,y) £ R+ x R™ | x > 0,e T x = l,y > 0} and So = {(x,y) £ S | ||t/|| < 1}. 

It is clear that So is a compact convex subset of S. Define the function F : S x S —> R 
by 

F(x, y, z, w) ={-Bx m - f(x,y)Ax m ~ 1 + (m - l)^diag (y)x [m ~ 2 \z) 

+ (i m - l)™~ 1 x lm ~ 1] - f(x,y)y [m ~ 1] ,w), (4.2) 

where 

_ m{m-l)^- 1 y T x^ rn - 1 ^ - Bx m 
iy x ^y) - Ax m + y T yl m ~ 1 ] 

Apparently, F(x, y, x, y) = 0 holds for any (x, y) £ S. Moreover, it can be seen that 
F(-,-;z,w) is lower-semicontinuous on S for any fixed (z,u>) £ S, and F(x,y,-,-) 
is concave on S for any fixed ( x,y ) £ S. With the given condition (4.1), we claim 
that 


fi := {(z,w) £ S | F{x,y-z,w) < 0, V (x,y) £ S 0 } 
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is compact. Otherwise, if 17 is not compact, then exists a sequence {( 2 ^, 11 /*^)} 
of 17 such that 

u/ fc ^)|| —> +00 as k —> + 00 . 

Since { 2 ^} is bounded, without loss of generality, we claim that ||«/*4|| —> + 00 . 
As a consequence, there exists io £ [n] such that w\ k ^ —¥ + 00 . By taking = 
= ei 0 £ So with ei 0 being the *o-th unit vector in R n , we have 


F(*«j,< fc >;2< fc \ «,<*>) = 6 k + K-io + 1 ~ m)(m - 1) 

a in...in i 1 


”*_ft ip ■■■ip ...(fc) 


where 

6 k = {-B{x^) m - 1 -f(x ik \y^ ) )A(x^ k) ) m - 1 +{m-l)^dia S (y ik) )(x ik) ) [m - 2] ,z ik) }. 
Clearly, the sequence {9 k } is bounded. It follows from the condition ( 6 ) that 

F(x^ k \y^-z^,w^)> 0 

for enough large k, which contradicts the fact that (z^ k \w^) £ 17. By Theorem 
6 in [13], there exists ( x,y ) £ S such that 

F(x,y-, z,w) <0, V (z,w) £ S. (4-3) 

Take w = 0 in (4.3), we know that, for any 2 £ D := {z £ R n | 2 > 0, e T z = 1}, 

F{x, y ; 2 ,0) = (—— fAx ml + (m — 1) ™ diag(i/)a^ m_2 l, 2 ) < 0, 

where / = f(x,y ), which implies 

Bx m ~ 1 + fAx m l — (m — l)”*diag (y)x^ m ~ 2 ^ > 0, (4-4) 

since D is a basis of R". Take again any w £ R”, it is clear that (x, y + w) £ S. 
Consequently, it holds that 

F(x,y;x,y + w) = F(x,y-,x,y) T ((m- - fy [m ~ 1] ,w) 

= (( m ~ - fy [m ~ 1] ,w) 

< 0 , 

where the second equality is due to the fact that F(x,y;x,y) = 0. Hence, 
fy [m ~ 1] - (m- l )™- 1 * 1 ™- 11 > 0 . 

Since x > 0 and x ^ 0, there exists io £ [n\ such that Xi 0 > 0. Accordingly, 

7-m—l \ / -1 \ —— 1-m—l . cy 

fVi 0 >{rn-l)™ x io >0, 


(4.5) 
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which implies / > 0. Denote I(y) = {i £ [n] | yi = 0}. It is clear that I(y ) is a 
proper subset of [n]. By (4.5), it is obvious that = 0 for any i € I{y)- So 

fvT -1 = (m~ (4.6) 

For any i G [n]\I(j/), taking w = te; with t 6 R, it follows from i/i > 0 that 
(i,f + »)eS for any real number t with enough small |t|. Recalling (4.3), we have 

F(x,y-,x,y + w) < 0, 

which implies that 

((m- 1 )^- 1 xr 1 -/r 1 )i<0 

for i G [n]\I(y) and any real number t with enough small |t|. We immediately 
obtain 

(m- = jy?~\ V i G [n]\I(v). (4.7) 

By (4.6) and (4.7), it holds that 

fy [m ~ 1] = (m - l)-- 1 * 1 ™- 11 , (4.8) 

or equivalently, 

/svij=( m -l)-sj. ( 4 . 9 ) 

Combining (4.4) and (4.8) leads to 

0 < Bx m ~ 1 + ]Ax m ~ x - f-^x [rn ~ 1] 

= f~^ {f^Bx™- 1 + ]^Ax m ~ x - x [m ~ 1] } , 

which implies 

A m Ax m ~ 1 + A Bx™- 1 - x [m ~ 1] > 0, (4.10) 

where A = /™- 1 . Now we verify that 

(x, A m Ax m ~ 1 + XBx m ~ 1 - x [m ~ 1] \ = 0. 

We only need to verify 

(xjAx™- 1 +Bx rn ~ 1 - =0, 

that is, 

n 

fAx m + Bx m - = 0. 

i=l 

Since F(x,y;x,y) = 0, that is, 


fAx m + Bx r " = m(m- l) m 11 - f^y™, 

i=l 
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we only need to further verify 

n n 

m{m - 1 )™~ 1 y T x [m ~ 1] - f^y™ - F"^ 1 = 0. (4.11) 

2=1 2=1 

Actually, the left hand of (4.11) amounts to 

n 

m(m- l)^- x y T x [m ~ 1] -fy T y lm ~ 1] ~ 

i=1 

n 

/ -i \ ^—1-T-fm—1] / \-1—1-T — \m —11 J-W \ ' -m 

= m(m — 1) m y x l J — ym — 1) m y x 1 ' — j m ~ 1 > Xj 

i=l 

/ -T —[m—11 J 11 

= (m — l) m y x —j m ~ 1 x x l 1 

= (m — l)™j/ T ^ m_ ^ — F m_1 i m — 1)"* /"*- 1 t/ T x^ m ^ 

= 0, 

where the first equality is due to (4.8), and the second equality comes from (4.9). 
Therefore, (A, x) is an m-degree Pareto-eigenpair of Q. □ 

Similarly, when we deal with the case of C := X, we can also establish the 
following result. 

Theorem 4.3 Let Q = (A, 13,1) £ T m ,n- Assume that A is strictly copositive. If A 
and B satisfy the following condition 

(m + au...i - 1 ){m - 1)™ _1 + bu...i >0, V i = 1, 2,... ,n, 

then Q has at least an m-degree Pareto-eigenpair. 

Proof Define the function h : R” x R" —> R by 

v ’ Ax m + y T y\. m ~d 

and the function G : S x S —> R by 

G(x,y; z,w) = (—Bx rn ~ 1 — h(x,y)Ax rn ~ 1 — (m—l)^dia.g(y)x^ rn ~ 2 \z) 
+({m-l)^- 1 x^-^-h{x,y)y^- 1 \w), 

where S is defined in the proof of Theorem 4.2. We can prove the assertion in a 
similar way that used in Theorem 4.2, and skip its details here. □ 

As a byproduct of Theorem 4.3, we immediately obtain the following existence 
result of the solution for QEiCP, which differs from the one presented in [38] . 

Corollary 4.1 Consider QEiCP corresponding to the special case of THDEiCP with 
m = 2. Let Q := ( A,B,I ) £ Mn- Assume that A is strictly copositive matrix. If A 
and B satisfy that an + bn + 1 > 0 for every i £ [n], then Q has at least one quadratic 
Pareto-eigenpair. 
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5 Numerical algorithm and experiments 

In this section, we first introduce an implementable splitting algorithm based upon 
the augmented Lagrangian method, which efficiently exploits the weakly coupled 
structure of the resulting optimization formulation of THDEiCP. Then, we conduct 
some computational results to show the reliability and convergence behavior of the 
proposed algorithm. 


5.1 The algorithm 

Note that model (3.1) can be recast as the standard minimization problem: 

min Bu m + #v T u [m - 1] 

s.t. Au m +lv m = l, (5.1) 

u > 0, v > 0, 

where d is a constant given by d := —m(m — 1)™ _1 . Here, we should notice that 
it is possible to employ the powerful semismooth and smoothing Newton methods 
to solve the model under consideration. However, we show below that a first- 
order structure-exploiting algorithm can be developed, which is much easier to be 
implemented than the second-order type methods. 

Taking a revisit on (5.1), we observe that (5.1) is an equality constrained opti¬ 
mization problem, and we know that the Augmented Lagrangian Method (ALM) 
[20,30] is a benchmark solver for this type model. Let C e R be the Lagrangian 
multiplier associated to the equality constraint. The augmented Lagrangian func¬ 
tion is given by 

3?{u,v,C) ■= Bu m +tiv T u [m - 1] -C{Au m +lv m - 1) + | (Au m +lv m - l ) 2 , (5.2) 

where /3 > 0 is the penalty parameter. Consequently, for a given £ R, the 
iterative scheme of ALM reads as follows: 

f (u ( - k+1 \v < ' k+1 ' ) ) = argmin j„S?(u, u,(^ fc) ) | m > 0, -u > 0 j ; (5.3a) 

| C (fe+!) = £(*) _ p (M(w (fc+1) ) m + l(v ik+1) ) m - l) . (5.3b) 

However, it seems not easy enough to implement such an algorithm due to the cou¬ 
pled structure and high nonlinearity emerging in the objective function and equal¬ 
ity constraint. To improve its implementability and numerical performance, the 
so-called Alternating Direction Method of Multipliers (ADMM) [17,19] was ju¬ 
diciously developed for separable convex minimizations by updating the variables 
in an alternating (Gauss-Seidel) order. In recent years, it is well documented that 
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ADMM has a surge of popularity in the areas such as signal/image processing, 
statistical learning, data mining, and so on. Here we just refer to [5, 12,18] for 
some surveys on ADMM. 

Following the spirit of ADMM, we split the first subproblem (5.3a) into two 
parts. For given we immediately have the following ADMM scheme: 



= argmin j„Sf(u, | u > 0 j ; 

(5.4a) 

u( fe+1 > 

= argmin ^J^{vS k+1 \v, C^) \ v > o| ; 

(5.4b) 

1 ) 

= C (fc) - P (A(u {k+1) ) m +l(v ( ' k+1) ) m - l) . 

(5.4c) 


It seems that such an algorithm exploits the weakly separable structure of model 
(5.1). However, it also fails to be easily implemented, because the first two sub¬ 
problems are not easy enough to have closed-form solutions. Indeed, we can clearly 
observe that both subproblems (5.4a) and (5.4b) have very simple convex sets as 
their constraints, thereby making the projections onto these sets very easy. Hence, 
it would greatly simplify (5.4) if both subproblems could reduce to the computa¬ 
tion of projections. Below, we consider the linearized version of (5.4) so that each 
subproblem has closed-form representation. Since (u, v, £) is nonconvex with 
respect to u and v in general cases, for the purpose of making both subproblems 
well-posed, we attach two proximal terms ^ ||« — iS k ^ || 2 and ^\\v — \\ 2 to (5.4a) 

and (5.4b), respectively. Here 71 and 72 are two positive constants. More specifi¬ 
cally, linearizing the nonlinear parts of Jz?(u, v , £) (see gradients in (3.2)), we derive 
a linearized ADMM as follows: 


u (fc+1) = 17 k - 

r (fc) * (fc) i 

u K ’ - 



7i J 


v {k+1) = 17 R? 

(k) ^( U < ~ k+1) ) [m - 1] +r( fc )' 



72 


C (fe+ i) = C (fc) - P [A{u {k+1) ) m +l(v {k+1) ) m - l) , 


(5.5a) 

(5.5b) 

(5.5c) 


where 17ij [■] represents the projection onto R"; 

:=mB(«®) ra 1 + i9(m — l)diag(i/ fc ^)(u^)[ m 2 ^+/3mq^A ( u ( fc )) m—1 
with qW := M(w (fc) ) m + l(v W ) m - 1 - and 


- (fe) :=/3 m ^4(u (fe+ 1 ) ) m +X( 


(k)\m 1 C ( } \ ~r[ (k)\m— 1 

- 1 - — 5 - X(ih >) 


Obviously, the linearized version (5.5) is more implementable than (5.4) due to 
the pretty simple iterative scheme. To the best of our knowledge, there is no 
convergence result of such a linearized ADMM for solving the underlying nonconvex 
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model. Therefore, it seems that our method (5.5) goes beyond the theoretical 
guarantees of the traditional ADMM. However, we will illustrate that our method 
(5.5) indeed is numerically convergent for model (5.1) in many cases. 


5.2 Numerical experiments 


We have shown theoretically that THDEiCP (1.3) is solvable when K := R” in 
Section 4 and introduce an implementable splitting method in Section 5.1. Now, 
we turn our attention to verifying our theoretical results and convergence behav¬ 
ior of the proposed algorithm (5.5) through preliminary computational results. We 
implement our algorithm by Matlab R2012b and conduct the numerical exper¬ 
iments on a Lenovo notebook with Intel(R) Core(TM) i5-5200U CPU@2.20GHz 
and 4GB RAM running on Windows 7 Home Premium operating system. 

Notice that model (5.1) is also available for matrix cases, and it is a new 
formulation in the QEiCP literature. Thus, we here also test a matrix scenario for 
the purpose of showing the efficiency of our proposed algorithm (5.5) in solving 
QEiCP. In the following experiments, we test three synthetic examples and only 
list the details of A (or A) and B (or B) in the coming examples. For the random 
data, we generate them by the Matlab script ‘rand’. 


Example 5.1 This example considers a special case of THDEiCP, that is, QEiCP, 
whose matrices A and B are uniformly distributed in (0,1) and given by 


A = 


t 0.2296 0.6870 0.7421 0.8943 \ 
0.6870 0.9403 0.1194 0.5919 
0.7421 0.1194 0.9325 0.7779 
v 0.8943 0.5919 0.7779 0.3290 j 


B = 


( 0.2235 0.3014 0.7879 0.5394 \ 
0.3014 0.4026 0.5329 0.5453 
0.7879 0.5329 0.8272 0.5375 
^ 0.5394 0.5453 0.5375 0.5994y 


Example 5.2 We consider the case where A and B are two 3-rd order 4-dimensional 
tensors; A is strictly copositive, but not nonnegative, and B is a randomly gener¬ 
ated tensor, whose entries are uniformly distributed in (1,2), that is, 


M(:,:,l) = 


/ 2 2 4/3 4/3 \ 

2 4/3 2/3 4/3 
4/3 2/3 8/3 0 


6 (:,:, 1 ) = 


( 1.6557 1.3572 1.7523 1.6055 \ 
1.3572 1.7577 1.4572 1.2192 
1.7523 1.4572 1.7060 1.0645 


A(:,:, 2) = 


y4/3 4/3 

0 

2 7 


\ 1.6055 1. 

2192 1. 

0645 1. 

8235 / 

( 2 

4/3 

2/3 

4/3 \ 



/1.6551 

1.5612 

1.4351 

1.6946 \ 

4/3 

12 

-2/3 

10/3 


:.2) = 

1.5612 

1.3404 

1.4202 

1.5916 

2/3 

-2/3 

16/3 

-2 

1.4351 

1.4202 

1.5060 

1.6231 

\ 4/3 

10/3 

-2 

14/3 ) 



y 1.6946 

1.5916 

1.6231 

1.1386 7 
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-40,:, 3) = 


4) = 


4/3 2/3 8/3 0 

2/3 -2/3 16/3 -2 
8/3 16/3 6 2/3 

0 -2 2/3 0 , 

^4/3 4/30 2 ^ 

4/3 10/3 -2 14/3 
0 -2 2/3 0 

2 14/3 0 4 j 


3) = 




1.9172 1.3331 1.6440 1.6613 
1.3331 1.5678 1.4275 1.2617 
1.6440 1.4275 1.9340 1.0709 
v 1.6613 1.2617 1.0709 1.3371, 

1.5383 1.5514 1.4477 1.3513 ^ 
1.5514 1.9619 1.4367 1.7875 
1.4477 1.4367 1.0844 1.4156 
1.3513 1.7875 1.4156 1.9106 , 


Example 5.3 We consider two 4-th order 3-dimensional symmetric tensors A and 
B, where A and B are randomly generated and uniformly distributed in (0,1). 
Specifically, A and B are taken as follows: 


0.6229 0.2644 0.3567 


0.3567 0.7367 0.1259 


-4(: fr :, 1,1) = 0.2644 0.0475 0.7367 ,B(:, :,!,!) = 0.4018 0.9957 0.0483 


0.6954 0.4018 0.1406 


0.1406 0.0483 0.0988 


0.7563 0.5878 0.5406 


0.5406 0.0715 0.3725 


-4(:,:, 1,2) = 0.5878 0.1379 0.0715 , 1,2)= 0.5351 0.2853 0.3071 


0.6730 0.5351 0.4473 


0.4473 0.3071 0.9665 


0.0657 0.4918 0.9312 


0.9312 0.9045 0.8711 


A {\,:, 1,3) = 0.4918 0.7788 0.9045 , 1,3)= 0.6433 0.8986 0.3427 


0.7585 0.6433 0.2306 


0.2306 0.3427 0.5390 


0.7563 0.5878 0.5406 


0.5406 0.0715 0.3725 


A(:, 2,1) = 0.5878 0.1379 0.0715 , i3(:,:, 2,1) = 0.5351 0.2853 0.3071 


0.6730 0.5351 0.4473 


0.4473 0.3071 0.9665 


0.7689 0.3941 0.6034 


0.6034 0.3465 0.4516 


A {\,:, 2,2) = 0.3941 0.3577 0.3465 , 2,2)= 0.3914 0.6822 0.5516 


0.3608 0.3914 0.5230 


0.5230 0.5516 0.7091 


0.8077 0.4910 0.2953 


0.2953 0.5556 0.9608 


A (:,:, 2,3) = 0.4910 0.5054 0.5556 2,3) = 0.2043 0.7282 0.7400 


0.4632 0.2043 0.2823 


0.2823 0.7400 0.9369 


0.0657 0.4918 0.9312 


0.9312 0.9045 0.8711 


A {\,:, 3,1) = 0.4918 0.7788 0.9045 , 3,1)= 0.6433 0.8986 0.3427 


0.7585 0.6433 0.2306 


0.2306 0.3427 0.5390 


0.8077 0.4910 0.2953 


0.2953 0.5556 0.9608 


A (:,:, 3,2) = 0.4910 0.5054 0.5556 , B(:,:,3,2)= 0.2043 0.7282 0.7400 


0.4632 0.2043 0.2823 


0.2823 0.7400 0.9369 
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/0.7581 0.7205 0.9044 \ /0.8200 0.5914 0.4983 \ 

A(:, :, 3,3) = 0.7205 0.0782 0.7240 , B(:,:,3,3)= 0.5914 0.0762 0.2854 . 

\ 0.9044 0.7240 0.3492/ \ 0.4983 0.2854 0.1266/ 

Before our experiments, we first introduce a reasonable stopping rule for the 
proposed method (5.5). Without loss of generality, we can use 

RelErr := max j||u (fc+1) — tt (fe) ||, \\v (k+1) - v (k) ||, |V (fe) |}<Tol (5.6) 

as a termination criterion to pursue an approximate solution with a preset toler¬ 
ance ‘Tol’, where \V^\ := |A(u (fc+ 1 ) ) m +X(v (fc+ 1 ) ) m -l| measures the violation of 
the underlying equality constraint. For the parameters involved in our algorithm, 
we throughout take f3 = 1, in addition to setting the starting points u ^ and 
v as randomly generated vectors and </ 0) = 0. For the other two parameters, 
we choose 71 = 200 and 72 = 10 for Example 5.1, 71 = 1000 and 72 = 50 for 
Examples 5.2-5.3. The tolerance ‘Tol’ in (5.6) is taken as Tol = 10 ~ 6 for all tests. 

As we have mentioned in introduction, tensor-related polynomial optimization 
problems suffers from high nonlinearity. From a theoretical point of view, lineariza¬ 
tion may destroy structural properties of the underlying functions, thereby result¬ 
ing in nonideal approximations so that the algorithm is not necessarily convergent 
for some cases. To investigate the performance of such a linearization, in Fig. 1, 
we plot evolutions of the relative error (‘RelErr’ defined by (5.6)) and the ob¬ 
jective value (‘Obj.’) of (5.1) [i.e., — , v^)] with respect to the number of 

iterations, respectively, and the ability of finding ideal solutions of THDEiCP. 



0 200 400 600 800 1000 1200 1400 1600 

number of iterations 



Exam. 5.1 Exam. 5.2 Exam. 5.3 


Fig. 1 Performance of the proposed algorithm. The left plot corresponds to the evolutions of 
‘RelErr’ defined by (5.6) and the objective value of (5.1) [i.e., with respect 

to the number of iterations, respectively. The right one shows the ability and reliability of the 
algorithm by testing 100 randomly generated starting points. 


From the left plots of Fig. 1, we can see that our linearized ADMM is convergent 
very fast with a random starting point. More importantly, it clearly shows that 
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the high nonlinearity leads to severely oscillating property in terms of the relative 
error. Actually, our computational experiences tell us that an inappropriate initial 
point far away a local solution may lead to divergence (see the right plot in Fig. 
1 ), which also implies that designing an implementable and stable algorithm for 
THDEiCP is a challengeable task. To further verify the ability and reliability 
of our algorithm, we randomly generate 100 different starting points such that 
u (0) = v ( 0 ) an( j en tries are uniformly distributed in (0,1). As we have proved 
in Theorem 3.1, a solution u^ of THDEiCP must be a nonzero vector. In our 
experiments, we observe that the algorithm terminates at a zero point in some 
cases. Accordingly, we record all results of the 100 tests and divide them into three 
groups: the first group corresponds to the divergent cases, which means the number 
of iterations exceeds the preset maximum iterations 20000 ; the second group refers 
to convergent cases but failed to find a nonzero solution of THDEiCP; the last 
group contains the cases of successfully finding a nonzero solution of THDEiCP. 
The corresponding rate of each group is graphically shown by the right plot in Fig. 
1 , which empirically exhibits the ability and reliability of our proposed algorithm. 
Indeed, we did a lot of experiments on QEiCP, and interestingly, all numerical 
results shows that the proposed linearized ADMM is always convergent for QEiCP. 
Thus, such an algorithm further enriches the solvers tailored for QEiCP. This 
also leaves us an open problem of whether there exists provable global or local 
convergence for the linearized ADMM on solving the nonconvex model (5.1). 

In Table 1, we list several groups of numerical results including starting points 
m*- 0 ) , eigenvalue (EigVal A), eigenvector (*), dual variable g, number of iterations 
(‘Iter.’) and computing time in seconds (Time), where the dual variable is defined 
by 

o •= A m Ax m ~ 1 + \Bx rn ~ 1 - Tx m ~ 1 

which together with the eigenvector x satisfies x 7 g = 0 . 

It can be easily seen from the data in Table 1 that our proposed algorithm is fast 
and reliable for solving the model under consideration, even it is not necessarily 
convergent for some cases. We conducted many simulations on random data and 
observed that an appropriate initial point and the parameters, especially the two 71 
and 72 are very important for convergence. Therefore, we will pay our attention 
on the study of convergence results of the iterative scheme (5.5) in the future. 
Additionally, we will provide some practical suggestions on choices of 71 and 72 . 


6 Conclusions 

This paper considers a unified model of THDEiCP including TGEiCP and QE¬ 
iCP as its special cases. As the first work on finding higher-degree cone eigen- 
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values of tensors, we analyze some corresponding topological properties including 
closedness, boundedness, and upper-semicontinuity of the m-degree /('-spectrum 
of Q (i.e., cr(Q,K) defined by (1.4)), and the number of m-degree Pareto- and K- 
eigenvalues of THDEiCP. For the special case where the underlying tensors A and 
B are symmetric, C = —X, and K := R™, we present a weakly coupled optimiza¬ 
tion formulation for THDEiCP, which is also a new formulation for QEiCP in the 
literature. Moreover, such a formulation could bring some numerical benefits for 
algorithmic design, for instance, an implementable linearized ADMM is developed 
in this paper. Theoretically, we establish results concerning existence of solutions 
of THDEiCP with general (symmetric and asymmetric) A and B when C := ±Z 
and K := M". In the future, we will study along this line, but with general ten¬ 
sor C not being a unit tensor (i.e., C ^ ±1). Of course, the convergence analysis 
of the proposed algorithm and designing algorithms for THDEiCP, especially for 
asymmetric cases, are also our future concerns. 
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